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ABSTRACT 

The growth of a supermassive black hole (BH) is determined by how much gas the host 
galaxy is able to feed it, which in turn is controlled by the cosmic environment, through 
galaxy mergers and accretion of cosmic flows that time how galaxies obtain their gas, 
but also by internal processes in the galaxy, such as star formation and feedback from 
stars and the BH itself. In this paper, we study the growth of a 10 12 M 0 halo at 
z = 2, which is the progenitor of a group of galaxies at z = 0, and of its central BH 
by means of a high-resolution zoomed cosmological simulation, the Seth simulation. 
We study the evolution of the BH driven by the accretion of cold gas in the galaxy, 
and explore the efficiency of the feedback from supernovae (SNe). For a relatively 
inefficient energy input from SNe, the BH grows at the Eddington rate from early 
times, and reaches self-regulation once it is massive enough. We find that at early 
cosmic times z > 3.5, efficient feedback from SNe forbids the formation of a settled 
disc as well as the accumulation of dense cold gas in the vicinity of the BH and starves 
the central compact object. As the galaxy and its halo accumulate mass, they become 
able to confine the nuclear inflows provided by major mergers and the BH grows at a 
sustained near-to-Eddington accretion rate. We argue that this mechanism should be 
ubiquitous amongst low-mass galaxies, corresponding to galaxies with a stellar mass 
below < 10 9 M 0 in our simulations. 

Key words: galaxies: formation — galaxies: evolution — galaxies: active — methods: 
numerical 


1 INTRODUCTION 

Supermassive black holes (BH) are common compact ob¬ 
jects observed in the centre of galaxies. Their are sus¬ 
pected to grow along with their host galaxy as observa¬ 
tions suggest a strong scaling between BH masses and galaxy 
properties (e.g. Magorrian et al. 1998; Tremaine et al. 2002; 
Haring &; Rix 2004; Kormendy & Ho 2013). Gas accretion 
drives most of the mass growth of black holes at high- 
redshift, and this accretion mechanism returns a fraction of 
the rest-mass accreted energy into effective feedback for the 
host-galaxy, which in turns can explain the observed scal¬ 
ing relation between black holes and galaxies (Silk & Rees 
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1998; King 2003; Wyithe & Loeb 2003). Such a behav¬ 
ior has been successfully implemented in modern cosmo¬ 
logical simulations and the impact of active galactic nu¬ 
clei (AGN) feedback on the galaxy mass content, the cir- 
cumgalactic gas properties and the halo mass distribution 
can be dramatic (Di Matteo et al. 2005; Croton et al. 2006; 
Bower et al. 2006; Sijacki et al. 2007; Booth & Schaye 2009; 
Dubois et al. 2010, 2012; Teyssier et al. 2011). 

AGN feedback is supposed to have a strong impact 
on the stellar and gaseous content of the most massive 
galaxies, helping to solve the so-called ‘over-cooling’ prob¬ 
lem, i.e. the expected fast cooling times and consequently 
high star formation rates predicted for galaxies, at odd 
with observational results. This was for instance high¬ 
lighted by comparing theoretical models of galaxy evolution 
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to the observed galaxy mass function (Croton et al. 2006; 
Bower et al. 2006). Theoretical models that do not include 
any form of negative feedback from stars or AGN over¬ 
predict the mass function at both the low- and high-mass 
end. Including stellar feedback, namely through supernova 
(SN) explosions, mitigates the problem at the low-mass end, 
but the SN energy input is not sufficient to eject gas from 
galaxies with virial halo velocity above a critical value of 
~ 100 kms -1 (Dekel &; Silk 1986), nor affect the cold stream 
inflow in high redshift galaxies (Powell et al. 2011). The en¬ 
ergetics involved in growing BHs are in principle compara¬ 
ble to, if not higher than, the binding energy of the host 
galaxy (e.g. Silk V Rees 1998; Fabian 2012). If a fraction 
of the rest-mass accreted energy can be tapped, it can pro¬ 
foundly affect the gas distribution in and around the galaxy. 
AGN-driven outflows can push the gas in the outskirts of the 
halo, even impacting the cold streams (Dubois et al. 2013; 
Costa et al. 2014). It will change the gas properties by in¬ 
creasing its temperature and decreasing its density, thereby, 
delaying the gas inflows onto galaxies and will impact the 
star formation rate of the host galaxy. 

The masses of BHs at the center of galaxies are 
observed to scale linearly with their host bulge with a 
BH to bulge mass ratio of Mbh/M b ~ 10 -3 . Assum¬ 
ing that the amount of gas mass impacted by the cen¬ 
tral AGN wind scales linearly with the bulge mass, the 
wind velocity u w oc LAGN/M gas oc Mbh/Mo is constant 
for BHs accreting at Eddington (Lagn oc Mbh). There¬ 
fore, the AGN wind velocity should be independent of the 
galaxy mass and AGN feedback have a non-negligible im¬ 
pact on high-mass galaxies as well as low-mass galaxies 
provided that the accretion of gas on the BH is unar¬ 
rested. Measurements of the BH mass correlation with the 
pseudobulge mass show a larger scatter and fall below 
the correlation between BH mass and classical bulge re¬ 
lation (Hu 2008; Kormendy et al. 2011; Sani et al. 2011). 
A similar deficit in BH mass is also observed for low-mass 
BHs (Greene et al. 2010; Jiang et al. 2011; Graham V Scott 
2015). SN feedback that is an important driver of gas- 
mass loss in dwarf galaxies (e.g. Springel V Hernquist 2003; 
Dubois V Teyssier 2008; Dalla Vecchia & Schaye 2008, 
2012; Hopkins et al. 2014) also turns dark matter cusps 
into cores (Gnedin V Zhao 2002; Mashchenko et al. 2006; 
Pontzen V Governato 2012; Teyssier et al. 2013), and re¬ 
duces the depth of the gravitational well. Therefore, the 
SN feedback activity in low-mass galaxies can probably in¬ 
hibit the activity of the central BH by depleting the cen¬ 
tral gas reservoir, and limiting the impact of AGN feed¬ 
back. Booth & Schaye (2013) performed a numerical study 
of the interplay between SN and AGN feedback in a large- 
volume cosmological simulation (see also Sijacki et al. 2007; 
Crain et al. 2015; and Newton & Kay 2013 for detailed iso¬ 
lated and merger simulations. They found that SN and AGN 
feedback weaken each other’s strength when jointly included 
in their models, compared to including only one of them at 
a time. However, they do not discuss whether SN feedback 
is able to impact the growth of BHs. 

This paper is part of a series of three papers investi¬ 
gating the connection between SN feedback and BH growth 
(paper I, this paper), SN-driven gas turbulence and BH spin 
evolution (Dubois et al. 2014, paper II), and the impact of 
cosmic galaxy evolution on the spin of central supermassive 


BHs (Dubois et al. 2014, paper III). The main purpose of 
this paper is to investigate the mass evolution of the central 
BH of a zoomed cosmological halo at high redshift, called the 
Seth 1 suite of simulations, that have enough resolution to 
start resolving the substructures of the interstellar medium 
(ISM) over several Gyrs. We test two implementations of SN 
feedback with the kinetic modelling of Dubois V Teyssier 
(2008), and the prescription from Teyssier et al. (2013) in 
order to probe how much the strength of SN activity can 
alter the growth of the central BH. 

In Section 2 we introduce the initial conditions and the 
numerical models for the physics of galaxy formation, the 
BH mass growth and its associated AGN feedback. In Sec¬ 
tion 3 we detail our results on BH mass evolution in the Seth 
simulations and how SN feedback influences BH growth. Fi¬ 
nally, in Section 4 we summarize our results. 


2 NUMERICAL SET-UP 

The simulations are run with the adaptive mesh refine¬ 
ment code RAMSES (Teyssier 2002). The evolution of the gas 
is followed using a second-order unsplit Godunov scheme 
for the Euler equations. The Harten-Lax-Van Leer contact 
(HLLC, Toro et al. 1994) Riemann solver with a MinMod 
total variation diminishing scheme to reconstruct the inter¬ 
polated variables from their cell-centered values is used to 
compute fluxes at cell interfaces. Collisionless particles (DM, 
star and BH particles) are evolved using a particle-mesh 
solver with a cloud-in-cell interpolation. 

2.1 The Seth suite of simulations 

We assume a ACDM cosmology with total matter density 

= 0.3, baryon density D& = 0.045, dark energy den¬ 
sity Da = 0.7, amplitude of the matter power spectrum 
cr 8 = 0.8285, n s = 0.9635 spectral index and Hubble con¬ 
stant Ho — 68.14kms -1 Mpc -1 consistent with the Planck 
data (Planck Collaboration et al. 2014) for the initial con¬ 
ditions produced with MUSIC (Hahn & Abel 2011). The box 
size of our simulation is Lbox — 50h _1 Mpc, with a coarse 
grid of 256 3 dark matter (DM) particles corresponding to 
a DM mass resolution of M res ,coarse = 8x 10 8 M©. A high- 
resolution region is defined around a halo of M v i r = 10 12 M 0 
at z = 2 that contains only high-resolution DM parti¬ 
cles within 2 r v ir (tv ir = 100 kpc) with mass M res ,hi g h = 
2 x 10 5 M 0 . The halo is a progenitor of a group of galaxies 
whose mass is M v i r = 5 x 10 12 M 0 at z = 0. 

The mesh is refined up to Ax = 8.7 pc (maximum level 
equals / max — 21 at redshift z = 3) or Ax — 34.8 pc (/max — 
19 at redshift z = 3) using a quasi-Lagrangian strategy: 
when more than 8 DM particles lie in a cell, or if the baryon 
density is larger than 8 times the initial DM resolution. The 
minimum cell size is kept roughly constant in proper size 
with redshift, i.e. an additional level of refinement is added 
every a exp — n x 0.1 (where n — 1,2,4, 8 and a exp is the 
expansion factor of the universe) up to level 21 (or 19) at 
a exp — 0.33 (z = 2). 

1 This name refers to the Egyptian god of storms, disorder and 
violence. 
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Figure 1. Projected gas densities of the Seth simulation (SHDA run) for the zoom region (top panel), the central halo (left panel), and 
the central galaxy (right panel) at 2 : = 3.25. The solid circle indicates the virial radius r v i r of the dark matter halo, and the dashed circle 
is for 0.1r v i r . We can see several cosmic cold filaments that connect to the central halo with a a complex gas distribution once they reach 
the central galaxy. The galaxy clearly exhibits a rotating disc-like structure with dense star-forming clumps of gas. 


2.2 Physics of galaxy formation 

Gas is allowed to cool by H and He cooling with a con¬ 
tribution from metals using a Sutherland & Dopita (1993) 
model for temperatures above T 0 = 10 3 K, which is the 
minimum temperature of the gas allowed through radia¬ 
tive losses. Heating from a uniform UV background takes 
place below redshift z re i on = 10 following Haardt & Madau 
(1996). Metallicity is modelled as a passive variable for the 
gas advected with the flow (whose composition is assumed 
to be solar) and is altered by the injection of gas ejecta dur¬ 
ing SN explosions and stellar mass losses. We assume a zero 


initial metallicity. The gas follows an equation of state for 
ideal monoatomic gas with adiabatic index 7 = 5/3. 

The star formation process is modelled with a Schmidt 
law: p* = e*p/ts , where p* is the star formation rate den¬ 
sity, e* the constant star formation efficiency, and t& the 
local free-fall time of the gas. We choose a low star for¬ 
mation efficiency e* & 0.02 consistent with observations of 
giant molecular clouds (Krumholz & Tan 2007) and surface 
density relations of galaxies (Kennicutt 1998). Star forma¬ 
tion is allowed in regions where the gas density exceeds 
no = 250Hem -3 for Ax = 8.7pc (no — 15 Hem -3 for 
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Table 1. List of the Seth simulation runs. First column indicates 
the simulation name. Second column indicates the minimum cell 
size. Third column indicates which SN feedback model is used. 
Fourth column indicates wether AGN feedback is included. Fifth 
column indicates the final redshift of the simulation. See the text 
for details. 


Name 

Ax (pc) 

SN 

AGN 

^end 

SL 

34.8 

no 

no 

2 

SLD 

34.8 

delayed cooling 

no 

2 

SLDA 

34.8 

delayed cooling 

yes 

2 

SHDA 

8.7 

delayed cooling 

yes 

2.4 

SLK 

34.8 

kinetic 

no 

2 

SLKA 

34.8 

kinetic 

yes 

2 

SHKA 

8.7 

kinetic 

yes 

3.5 


Ax = 34.8 pc) and where the temperature is not larger 
than T 0 (corrected for the contribution of the polytropic 
equation of state, see below). The stellar mass resolution is 
77i s ,res = noAx 3 ~ 5 x 10 3 M© for the Ax = 8.7 pc resolu¬ 
tion runs (m s ,res — 2 x 10 4 M© for the Ax = 34.8 pc). The 
gas pressure is artificially enhanced above p > po assuming 
a polytropic equation of state T — To(p/po) K '~ 1 with poly¬ 
tropic index k = 2 to avoid excessive gas fragmentation. The 
value of no for both the density threshold of star formation 
and the polytropic equation of state is chosen so that the 
Jeans length is always resolved with 4 cells. 

Feedback from stars is taken into account assuming a 
Salpeter initial mass function with rj sn = 0.1 of the mass 
fraction of stars ending up into a type II supernova and 
releasing esN = 10 5 O ergMQ 1 . The numerical implementa¬ 
tion is done in two flavors, (i) kinetic : We use the imple¬ 
mentation from Dubois & Teyssier (2008) based on a Sedov- 
blast wave solutions: mass, momentum and energy are de¬ 
posited in the surrounding gas cells so that it mimics the 
Sedov solution with a mass loading of f w — m e j/m S;re s = 1, 
where m e j is the gas mass of the ejecta taken from the cen¬ 
tral cell and redistributed to all the cells within the bub¬ 
ble of size rsN = 1.5Ax (see Dubois & Teyssier 2008 for 
more details), (i) delayed cooling : An alternative modeling 
is that of Teyssier et al. (2013) which locally delays the gas 
cooling after an explosion (in the spirit of Stinson et al. 
2006), where energy is released both in the classical en¬ 
ergy component and in a tracer component, which is de¬ 
noted ‘non-thermal’ and is passively advected with the flow. 
The passive non-thermal energy decays on a typical time- 
scale of tdiss = 0.8 Myr for Ax = 8.7pc (tdiss = 2 Myr 
for Ax = 34.8pc), and the gas cooling is not allowed until 
the velocity dispersion associated to the non-thermal energy 
component is larger than ont = 50kms _1 . These values of 
tdiss and ctnt are different from the initial choice adopted 
in Teyssier et al. (2013) who set them up to the typical 
molecular cloud lifetimes (Williams & McKee 1997). Here, 
we prefer to adopt a resolution-dependent approach, so that 
the explosion is able to propagate over a few cell resolu¬ 
tion elements. The choice of ctnt and tdiss are discussed in 
Appendix A. 

2.3 Model for BH growth and AGN feedback 

We use the AGN feedback model described in Dubois et al. 
(2012). BHs are created at loci where both the gas den¬ 


sity and stellar density is larger than the density threshold 
for star formation po with an initial seed mass of M see d = 
min (0.75p gas Ax 3 , 10 5 M©) (i.e. between [1.5 x 10 4 ,10 5 ] M© 
for Ax = 34.8 pc and [3.7 x 10 3 ,10 5 ] M© for Ax — 8.7pc). 
Thus, the initial seed mass can be smaller than 10 5 M© 
(BH mass through direct collapse, Begelman et al. 2006; 
Spaans & Silk 2006) when gas has not accumulated enough 
to contain this amount of gas mass but has formed enough 
stars to reach the stellar density criterion. The accretion 
rate onto BHs follows the Bondi-Hoyle-Lyttleton (Bondi 
1952) rate Mbh = 4:7vaG 2 M^ u p/(c 2 s + u 2 ) 3 / 2 , where G is 
the gravitational constant, Mbh is the BH mass, p is the 
average gas density, c s is the average sound speed, u is 
the average gas velocity relative to the BH velocity, and 
a is a dimensionless boost factor with a — ( p/po) 2 when 
p > po and a = 1 otherwise (Booth & Schaye 2009) in or¬ 
der to account for our inability to capture the colder and 
higher density regions of the ISM. Averaged quantities are 
summed over 4 cells in radius and kernel-weighted with a 
Bondi radius-dependancy (see Dubois et al. 2012 for the 
exact form of the averaging procedure). The effective ac¬ 
cretion rate onto BHs is capped at the Eddington accre¬ 
tion rate: M E dd = 47rGMBum p /(e r aTc), where ctt is the 
Thompson cross-section, c is the speed of light, m p is the 
proton mass, and e r is the radiative efficiency, assumed to 
be proportional to the spin of the BH with e r = 1 — Ti sco = 
1 — y/l — 2/(3n S co). Spin magnitude and direction are fol¬ 
lowed self-consistently, on the fly, in this simulation. We refer 
the reader to paper II (see also Dubois et al. 2014) for details 
and we just summarize here the main assumptions. The spin 
magnitude of BHs is modified by accretion of gas through 
the expression derived by Bardeen (1970). If the angular mo¬ 
mentum of the accreted gas is misaligned with the direction 
of the BH spin, Lense-Thirring effect causes the innermost 
parts of the disc to rotate within the equatorial plane of the 
BH and a warped disc is created that eventually completely 
aligns or anti-aligns with the BH spin. The case for anti¬ 
alignment of the BH with the disc requires that (King et al. 
2005) cos# < — 0.5Jd/JBH, where J a and Jbh are the angu¬ 
lar momenta of disc and BH respectively, and 6 is the angle 
between Jbh and Jd- To model the unresolved accretion 
disc we adopt the Shakura & Sunyaev (1973) thin accretion 
disc solution, and define Ja to be the angular momentum cal¬ 
culated at smaller between the warp radius, i.e., the radius 
within which gas can effectively transfer its angular momen¬ 
tum to the BH within a viscous timescale, or the self-gravity 
radius, i.e. the radius beyond which the disc becomes prone 
to fragmentation. Finally, we update BH spins after BH-BH 
coalescences using the analytical fit of Rezzolla et al. (2008) 
derived from relativistic numerical simulations of BH bina¬ 
ries. 

In order to avoid spurious oscillations of the BH in the 
gravitational potential well due to external perturbations 
and finite resolution effects, we introduce a drag force that 
mimics the dynamical friction exerted by the gas onto a 
massive particle. This dynamical friction is proportional to 
Tdf = /gas 47 r ap(GMBH /c s ) 2 , where / gas is a fudge factor 
whose value is between 0 and 2 and is a function of the Mach 
number M = u/c s < 1 (Ostriker 1999; Chapon et al. 2013), 
and where we introduce the boost factor a for the same rea¬ 
sons than stated above. A complementary approach, that we 
do not model here since gas dynamical friction is sufficient, 
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is to explicitly model the unresolved dynamical friction ex¬ 
erted by DM particles onto BHs (Tremmel et al. 2015). 


The AGN feedback is a combination of two differ¬ 
ent modes, the so-called radio mode operating when x — 
Meti/MEdd < 0.01 and the quasar mode active otherwise. 
The quasar mode corresponds to an isotropic injection of 
thermal energy into the gas within a sphere of radius Ax, at 
an energy deposition rate: ^agn — CfCr-MsHC 2 , where ef = 
0.15 for the quasar mode is a free parameter chosen to re¬ 
produce the M BH -M b , MsH-cTb, and BH density in our local 
Universe (see Dubois et al. 2012). At low accretion rates on 
the other hand, the radio mode deposits the AGN feedback 
energy into a bipolar outflow with a jet velocity of 10 4 km s -1 
into a cylinder with a cross-section of radius Ax and height 
2 Ax following Omma et al. (2004) (more details about the 
jet implementation are given in Dubois et al. 2010). The ef¬ 
ficiency of the radio mode is larger with ef = 1. Though this 
quasar/radio mode matches the observed radiatively effi¬ 
cient/inefficient activities of AGN feedback (Churazov et al. 
2005; Russell et al. 2013), it is possible that some accretion 
rate-driven jets becomes indistinguishable from large open¬ 
ing angle AGN-driven nuclear outflows after they have prop¬ 
agated into the multiphase gas of the galaxy (Wagner et al. 
2012 ). 


The list of the Seth simulation runs is listed in Ta¬ 
ble 1. All simulations with Ax = 34.8 pc resolution are 
run up to z = 2, while the simulations with Ax =. 8.7 pc 
resolution are stopped at more intermediate redshift after 
the BH of the central simulation galaxy has reached self- 
regulation. We perform a simulation with no feedback from 
SNe nor from AGN (and without BH formation) denoted SL 
(Ax = 34.8pc). We have two simulations with SN feedback, 
one version with kinetic implementation SLK and another 
one with delayed cooling SLD, which do not include BHs 
(Ax = 34.8pc). Finally, we run four simulations including 
AGN feedback, two with kinetic SN feedback SLKA and 
SHKA (respectively with Ax = 34.8pc and Ax = 8.7pc), 
and two with the delayed cooling SN feedback SLDA and 
SHDA. 


Fig. 1 shows three projections of the gas density for the 
SHDA run at z = 3.25, one at large scale where cosmic fila¬ 
ments are clearly visibles (within a large-scale cosmic wall), 
one at the halo scale where these cosmic filaments penetrate 
the central halo in the form of dense and cold gas mixing 
with the diffuse ambient medium, and another projection at 
galactic scale where we see the detailed structure of the gas 
within the galaxy in the form of dense star-forming clumps 
embedded into gas spiral arms. In those projections, we have 
indicated the virial radius r v n of the dark matter halo (as 
the radius at which the average density in the halo is 200 
times the average cosmic density). With zoom simulations, 
we can resolve at the same time the complex pattern of the 
large-scale inflows and outflows of gas into the halo and onto 
the galaxy, and the multiphase structure of the interstellar 
gas. 
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redshift 

14.0 8.5 6.2 4.9 4.1 3.5 3.1 2.7 2.4 2.2 2.0 



time (Gyr) 

Figure 2. Top: BH mass as a function of time (redshift) for SHDA 
(black), SLDA (grey), SHKA (red) and SLKA (orange). Bottom: 
BH mass accretion rate (solid lines) and Eddington mass accre¬ 
tion rate (dotted lines) for SLDA and SLKA simulations only and 
values are smoothed for sake of readability with a 5 Myr time win¬ 
dow. The small inset shows the exact BH mass accretion of the 
SLDA simulation between 2.1 <t< 2.2 Gyr. The horizontal lines 
indicate the SN and AGN regulation phases of the BH growth 
with a color that corresponds to the simulation it refers to. The 
simulations with kinetic SN feedback (SLKA and SHKA) show a 
rapid growth of the BH mass at early times and quickly regulate 
the BH mass growth. The simulations with delayed cooling SN 
feedback (SLDA and SHDA), instead inhibit the early BH growth 
and delay its Eddington-limited phase to later times (below red¬ 
shift z < 3.5). 

3 RESULTS 

3.1 Quenching early BH growth with efficient SN 
feedback 

The top panel of Fig. 2 shows the BH mass evolution for 
the four simulations including AGN feedback, and the bot¬ 
tom panel of Fig. 2 shows the BH mass accretion rate for 
SLDA and SLKA only. Both simulations with kinetic SN 
feedback, SLKA and SHKA, show a rapid initial growth of 
the BH mass at the Eddington rate between 6 < z < 9, 
with a plateau below z < 6 when the BH has reached self- 
regulation. Viceversa, both simulations with delayed cooling 
SN feedback, SLDA and SHDA, have an extended period of 
time, down to z ~ 3.5, where the BH grows quiescently, 
far from the Eddington limit. It is only below this redshift, 
that the BH growth takes off at approximately half the Ed¬ 
dington rate from 2.5 < z < 3.5 and then settles in a self- 
regulated state where the BH mass is comparable in runs 
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redshift 


9.7 8.5 7.5 6.8 6.2 5.7 5.3 4.9 



Figure 3. Mass of gas within 1 per cent of the halo virial radius as 
a function of redshift for SHDA (black), SLDA (grey solid), SLD 
(grey dashed), SHKA (red), SLKA (orange solid), SLK (orange 
dashed), and SL (blue). Only simulations with delayed cooling 
SN feedback exhibit rapid variations in the gas mass content at 
early times due to the efficient evacuation of the cold gas from 
the galaxy center. 


SLDA, SHDA and SLKA. Note that the evolution of BH 
mass is relatively independent of resolution for the same set 
of SN sub-grid modeling: Eddington-limited phases appear 
at the same epoch and their time extensions are similar, and 
BH masses at self-regulation are the same to within less than 
a factor of 2. 

At z ~ 3.5, a significant galaxy merger with stellar mass 
ratio of-1:3 drives a large inflow of gas towards the BH. The 
three simulations that have reached this redshift (SLKA, 
SHDA and SLDA) show a rapid increase of the BH mass. In 
particular, the simulations with delayed cooling SN feedback 
have an extended near-Eddington growth for the central BH 
(2.5 < 2: < 3.5), due to the absence of BH self-regulation 
prior to the merger. Thus, the BH grows until it reaches a 
mass sufficiently large to pump enough energy into the gas 
to unbind it from the galaxy gravitational potential well. 

As explained in detail in paper II, delayed cooling SN 
feedback causes cold gas clumps in the galaxy to be regu¬ 
larly destroyed by SN explosions, every —10 Myr. Periods of 
intense star formation are followed by phases where stellar 
feedback drives large-scale outflows and hot bubbles, and 
removes cold gas from the galaxy at the same time. During 
these phases, the BH accretion rate is well below the Edding¬ 


ton limit 2 , while it is at Eddington when the gas cools down 
again and re-accumulates around the BH. This behaviour is 
highlighted in the inset in Fig. 2. The net effect is to overall 
reduce the Eddington-limited growth phase by a factor of 2, 
as reflected by the time-smoothed curve (for sake of clarity) 
shown in the bottom panel of Fig. 2: the average BH mass 
accretion rate is a factor 2 below the Eddington accretion 
rate between 2.5 < z < 3.5 for the SLDA simulation (as well 
as for the SHDA, not represented here). 

In summary, the delayed cooling implementation of SN 
feedback quenches the initial growth of the BH, while the 
kinetic SN feedback allows for a rapid Eddington-limited 
growth at early times. The fundamental reason is that the 
two different approaches to SN feedback affect the den¬ 
sity and temperature of gas in the galaxy differently. Ki¬ 
netic SN feedback imparts only momentum to the gas 
(with lower magnitude than the expected analytical solu¬ 
tion, see Kimm et al. 2015). The gas quickly thermalizes if 
the conditions to create a shock are met, and such condi¬ 
tions are effectively met each time there is a SN explosion. 
This can be understood as follows. In the early phase, SN 
explosions are similar to Sedov solutions, therefore in a few 
time-steps, the solution will approach that of the adiabatic 
expansion phase, where 20 per cent of the energy is in kinetic 
form and 80 per cent in internal energy (e.g. Cox 1972). In 
dense star-forming regions, the temperature reached after a 
SN explosion is not high enough to avoid over-cooling of the 
internal component since the Sedov scale remains unresolved 
at our spatial resolution, and 80 per cent of the total energy 
of the explosion is lost, limiting the total amount of mo¬ 
mentum that would have been transferred by the hot pres¬ 
surised bubble by the end of the energy-conserving phase. 
Thus, if gas is not artificially prevented from cooling, dense 
gas can survive SN explosions, contrary to expectations. De¬ 
layed cooling SN feedback obviates this shortcoming by ex¬ 
plicitly forcing the gas to not cool in the dense regions of 
star-format ion, mimicking the effect of internal energy de¬ 
position (however, it does not guarantee that the correct 
amount of momentum is imparted to the gas, only that en¬ 
ergy is efficiently transferred to the gas). This is the reason 
why simulations with delayed cooling SN feedback show that 
dense star-forming gas is destroyed after the ignition of the 
first stars into SNe (see paper II). 

Fig. 3 clearly shows these different behaviours. The gas 
mass within 1 per cent of the virial radius (~ 100 pc at 
2: = 6) has a relatively smooth time-evolution with the ki¬ 
netic mode of SN feedback (SHKA, SLKA and SLK runs), 
comparable to the simulation without SNe (SL run). Simu¬ 
lations with delayed cooling SN feedback (SHDA, SLDA and 
SLD runs) instead show strong variations in the amount of 
gas mass relatively close to the BH. These variations are an 
intrinsic feature related to the modelling of SNe as they are 
detected independently of the presence of the AGN feedback 
and the choice of resolution. 

Fig. 4 highlights two specific times representative of the 

2 Note how the Eddington accretion rate exhibits negative vari¬ 
ations with time, while the BH mass can only increase. This ap¬ 
parent paradox is explained by variations in the spin parameter 
of the BH which in turn change the Eddington accretion rate by 
modifying the radiative efficiency e r (see paper II and paper III 
for a complete overview of the BH spin evolution). 


© 0000 RAS, MNRAS 000, 000-000 
















SN-regulated BH growth 


7 


SHKA 


edge—on 

SHKA 


edge—on 

/ 

/ 


\ 

/ 


. 

\ 

\ 






i 

\ 


/ 

\ 


/ 

1 

■ 

M ! z=5.80 

1 kpc 


1 z=5^5 


■ j 

1 1 


fenr i 




II Si 



SHKA 


face— on 

SHKA 


face—ort; 

■ 


. 



. 

■ 


f 

* 

1 

j « 

'Ssj? 

"\ 

i 

\ 


/ 

F 


/ 

1 kpc 


|i z=5.80 

1 kpc 

z=5.65 



Figure 4. Gas density projections of the central galaxy at z = 5.8 (first and third columns) and z = 5.65 (second and fourth columns) 
for the SHDA (four left-most panels) and SHKA (four right-most panels) simulations, seen edge-on (first row) or face-on (second row). 
The blue contours surround the regions of star formation with gas density above 250Hem -3 . The dashed circles indicate the 10 per cent 
virial radius of the halo. In the simulation with the delayed cooling SN feedback, dense star-forming clumps in the centre of the galaxy 
are regularly destroyed and reprocessed, while in the simulation with kinetic SN feedback the star-forming gas settles into a long-lived 
compact rotating disc. 


two phases characteristic of the delayed cooling SN feedback: 
intense star formation driven by cold gas followed by SN ac¬ 
tivity which in turns removes the cold dense gas. The four 
left-most panels of Fig. 4 show the projected gas density at a 
maximum (z — 5.80) and near a minimum (z — 5.65) of the 
gas mass within 1 per cent of the halo virial radius for the 
SHDA simulation (cf. Fig. 3, top panel). The dense gas in 
the galaxy appear edge-on with some cold clumps near the 
very center of the galaxy at z — 5.8. However, at z — 5.65, 
the dense gas has disappeared from the center of the image 
due to SN explosions that have disrupted this cold gas com¬ 
ponent. For comparison, the four right-most panels of Fig. 4 
show the same gas density projections for the SHKA simu¬ 
lation, where a thin disc of gas also appears nearly edge-on, 
and this disc is present in both time steps. Therefore, con¬ 
trary to the case with delayed cooling SN feedback (SHDA), 
the simulation with kinetic SN feedback (SHKA) does not 
show any sign of the destruction of the cold star-forming 
gas in the galaxy. Finally, the presence at all times of this 
cold gas reservoir near the BH for SHKA (and SLKA) is the 
reason for the rapid early growth of the BH, while in SHDA 
(and SLDA) disruption of dense gas is efficient enough to 
starve the BH. 

Fig. 5 shows the gas number density in the vicinity of 
the BH as a function of time, as well as the gas sound speed 
and its relative velocity with respect to the BH. During the 
early growth of the BH in the SLKA run (6 < z < 9), the gas 
density near the BH is above n > 10 3 cm -3 , the gas sound 
speed is ~ 50 kms -1 , and the gas relative velocity is u < 
10 kms -1 . Therefore, the gas is sufficiently dense, cold, and 
with very little turbulence (relative gas motion), to allow the 
Bondi accretion rate to be Eddington-limited, and the BH 


to grow to reach self-regulation in ~ 200 Myr. During the 
out-of-Eddington phases, gas density drops to lower values 
(n ~ 20 — 100 cm -3 ), the gas sound speed increases (c s ~ 
200 km s -1 ), as does the gas relative velocity (u ~ 100 — 
300 kms -1 ), all being consequences of the energy input from 
the AGN while the BH is at self-regulation. It is particularly 
striking that these quantities change abruptly (gas density 
increases and velocities decrease) with the new inflow of gas 
triggered by the merger at z ~ 3.5 (t ~ 2 Gyr). 

For the SLDA run, gas density, sound speed and gas rel¬ 
ative velocity values up to z = 3.5 are comparable to the val¬ 
ues obtained in the SLKA run during the BH self-regulation 
phase (i.e., after z ~ 6; t ~ lGyr). The main difference is 
that, in the SLDA run, the absence of Eddington-limited 
phase is due to the energy input from SNe that destroy the 
dense star-forming gas, keeping the density low, the temper¬ 
ature high and the interstellar medium turbulent. The gas 
number density in SLDA between 5 < z < 9 shows rapid 
variations with peaks at n ~ 500 cm -3 followed by drops to 
n ~ 50 — 100 cm -3 . It illustrates well the “breathing” of the 
galaxy when delayed cooling SNe feedback is employed (see 
also Fig. 3): dense gas accumulates in the galaxy (and near 
the BH) for At ~ 10 Myr until the first SNe explode and 
drive the cold gas in an expanding wind that falls back onto 
the galaxy in a dynamical time. At z cb 3.5, a merger trig¬ 
gers fresh gas near the BH that accumulates and leads to the 
near-Eddington phase: gas sound speed and gas relative ve¬ 
locity are still high (respectively ~300 and 100kms -1 ), but 
the gas density remains steadily at n ~ 500 cm -3 . Thus, 
even though velocities are on average large due to the star- 
burst activity induced by the merger, the gas density is suf¬ 
ficiently high to keep the BH growing near the Eddington 
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Figure 5. Gas number density (top), gas sound speed (middle), 
and gas relative velocity (bottom) near the BH as a function of 
time (redshift) for SLDA (grey) and SLKA (orange). 


rate. Note that these quantities have strong small-scale time 
variations (on a time-scale of 10 Myr) that cannot be rep¬ 
resented in Fig. 5. After the merger-driven triggering of the 
BH growth in the SLDA simulation, the BH ends up self¬ 
regulating at z ~ 2.5, and the gas density drops, while the 
gas sound speed and relative velocity increase to values com¬ 
parable to that of the SLKA simulation at the same epoch 
and for the same reason: its AGN activity is strong enough 
to remove the cold gas from its direct surroundings. 


3.2 Why SNe are able to quench BH growth in 
low-mass galaxies 

From our set of simulations, we have seen that the central 
BH growth is initially inhibited by the presence of SNe when 
they are associated with a delayed cooling prescription that 
maximises the gas motion. Now, we want to compare the 
theoretical expectation of the velocity of SN-driven winds 
to the escape velocity of the galactic bulge. Therefore, if 
SN-driven winds are efficient, they should drive the gas out 


redshift 

14.0 8.5 6.2 4.9 4.1 3.5 3.1 2.7 2.4 2.2 2.0 



time (Gyr) 

Figure 6. Escape velocity from the bulge (at r b = 100 pc from 
the galactic centre) as a function of time (redshift) for SL (blue), 
SLDA (grey solid), SLD (grey dashed), SLKA (orange solid) and 
SLK (orange dashed). The horizontal dot-dashed line corresponds 
to the velocity of SN-driven winds in the bulge with /blast — 1 
(corresponding to the delayed cooling SN model). Grey horizontal 
lines mark the SN- and AGN-dominated evolution in the delayed 
cooling case. The orange horizontal line marks AGN-dominated 
evolution in the kinetic feedback case. In SLDA and SLD, the 
end of the SN regulated growth of the BH corresponds to the 
time when the escape velocity of the bulge exceeds the velocity 
of SN-driven wind. In the SLKA and SLK runs, the BH grows as 
soon as it is created, and there is no SN-dominated phase, which 
indicates that the SN driven-winds for the kinetic model have a 
very low energy-conversion efficiency (/blast C 1) and, thus, an 
ejection velocity below the horizontal dot-dashed line. 


of the bulge (and star-forming clumps), until the bulge be¬ 
comes sufficiently massive to retain any outflows in its sphere 
of influence. 

The characteristic velocity of SN-driven winds in a star¬ 
forming clump is of the order of the Sedov velocity: 


usn 


1.2 


/ 


YTl new,s 


/blast ?7SNeSN 

ra g 


270 


/blast^SN / (m : 


0.1 


'r,s/ m s) kms -1 . 


0.1 


( 1 ) 


where esN — 10 50 ergM^ 1 is the SN type II specific energy, 
/blast ^ 1 is the fraction of SN energy input that is effectively 
transferred into the gas (in the case of efficient energy losses 
through gas cooling /blast can be much lower than 1), and 
ranew,s is the amount of new stars formed within the gas 
clump. ra ne w,s depends on the clump free-fall time (i.e. the 
clump gas density) and the time-scale before the ignition 
of the first SNe £sn = 10 Myr. Integrating the Schmidt law 
over time and volume, and assuming that the gas density is 
nearly constant over tsN, one obtains 


^new,s 

ra g 



( 2 ) 


where the gas free-fall time t& for the central gas clump 
of mass rsj 10 8 M® in the proto-bulge (see Fig. 3) and 
of ~ 100 pc radius 3 is ~ 1.6 Myr. Therefore, the ratio of 


3 This radius is obtained doing a bulge-disc decomposition on 
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Figure 7. Integrated mass profiles as a function of radius normalised to the halo virial radius at z = 4.71 (left panel) and z = 3.26 (right 
panel). Grey lines are for the SLDA simulation and orange lines are for the SLKA simulation. Solid lines are the total masses, dotted 
lines the DM masses, dashed lines are the stellar masses and dot-dashed lines are the gas masses. At early times z = 4.71, the simulation 
including SN feedback with delayed cooling (SLDA) has reduced the mass concentration within 10 per cent of the virial radius compared 
to the simulation with kinetic mode of SN feedback (SLKA). 


^new,s/% — 0.12. To summarize, the velocity of the blast 
propagating into the central gas clump (bulge) should be 
close to 270 kms -1 for the delayed cooling case, where 
/blast ~ 1 is ensured by construction. 

The escape velocity of a clump is 

_ / 2Gm c \ ^ 

^esc — \ I 5 ('ll 

V r c \ 

where m c \ and r c \ are respectively the clump total mass and 
clump radius. If we replace the clump radius r c 1 by the bulge 
radius r b — 100 pc, and the clump mass m c \ by the bulge 
mass Mb ~ 10 9 M© that are the representative values when 
the BH starts growing at the Eddington rate. The escape 
velocity is, therefore, u esc — 300 kms -1 and is comparable to 
the Sedov velocity produce by SN winds in the bulge. Thus, 
below Mb < 10 9 M©, and provided that /blast ~ 1, the gas 
cannot accumulate in the bulge and is regularly ejected from 
the bulge gravitational potential well (see the top panel of 
Fig. 3 for SHDA, SLDA and SLD), and above Mb ^ 10 9 M© 
the gas is locked into the bulge, in a fountain that ejects gas 
that quickly falls back, and provides the fuel to feed the 
central BH over a prolonged duration. 

Fig. 6 shows the escape velocity within 100 pc distance 
from the center of the galaxy. For sake of simplicity, we 
call that region the bulge here, but note, however, that the 
bulge is not always well defined, especially at early times 
when only short-lived star forming clumps are present in 

the stellar density profile of SHDA using a fitting procedure with 
two exponentials. The bulge scale radius is of ~ 25 pc. We use a 
radius of 100 pc in order to account for 90 per cent of the stellar 
mass content in the bulge in that case. 


the galaxy with delayed cooling SN feedback. At early times 
(z > 3.1), since SNe in the delayed cooling cases prevent the 
gas from accumulating too rapidly in the bulge, we observe 
that the escape velocity of the bulge increases more slowly 
than in the kinetic feedback cases. At z = 3.1, the escape 
velocity in the bulge in SLD and SLDA reaches 270 kms -1 . 
This value is the predicted SN-driven wind velocity in the 
bulge for /blast = 1 and marks the end of the SN-driven 
evolution of the bulge gas content. From this moment, the 
BH can grow close to the Eddington rate and settles into 
self-regulated growth ~ 400 Myr later at z = 2.6. Note 
that escape velocities in the bulge in both SLD and SLDA, 
respectively with and without AGN feedback, are similar 
during the SN-driven evolution since the AGN has as yet no 
impact on the surrounding gas as SN-driven winds prevent 
the BH from accreting gas at near-Eddington. There is a 
difference in escape velocity between SLDA and SHDA (i.e. 
when changing resolution) at late times. This is attributed to 
the more compact bulge in the high resolution run (~ 25 pc 
against ~ 60pc bulge scale radius), although the mass of 
the bulge remains identical (see section 3.3). 

For the kinetic cases, SLK and SLKA, the bulge en¬ 
dures a rapid mass growth already very early on due to 
the inefficiency of the coupling of the SN energy to the gas 
(/blast 1). Until the BH has reached its self-regulated 
phase in SLKA, the bulge mass can grow steadily, as it 
does in the case without any feedback from SNe nor AGN 
(SL run). Finally, the growth of the bulge escape veloc¬ 
ity (or equivalently the bulge mass) is suppressed once the 
BH is massive enough to have a significant impact on the 
surrounding gas, which is around z ~ 6 in SLKA. Inter¬ 
estingly, the bulge escape velocities of SLDA and SLKA 
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are similar u esc — 500kms _1 once both simulations have 
reached the BH self-regulated phase (z < 3), and their 
value is smaller compared to the runs without AGN feed¬ 
back Uesc — 800 km s -1 (SL, SLD and SLK). 

A more in-depth investigation of the structure of the 
gravitational potential well, and, thus, mass distribution is 
shown in Fig. 7, where we see how the gas mass accumulates 
differently in SLKA and SLDA at early epochs. In the SLKA 
run, the escape velocity peaks at ~ 400 km s -1 at z — 4.71 
within the galactic bulge (< 0.01 r v i r = 200 pc) due to the 
high concentration of stars in the central region of the galaxy 
(and with a low concentration of gas and DM). Vice versa, 
in the SLDA run at the same redshift, the peak of the escape 
velocity is 250 km s -1 at a distance ~ 0.05 r v i r = 1 kpc that 
roughly corresponds to the galaxy radius. The escape ve¬ 
locity within the galactic bulge reaches lower values around 
150 km s -1 a factor of two below the SN wind velocity. In 
this simulation, the dominant mass component in the bulge 
is the gas with a fraction above 50%. Though the amount of 
gas within the bulge is larger by a factor of 2 in SLDA than 
in SLKA, the amount of total mass is lower in SLDA because 
of the frequent expulsions of gas, that translate into a much 
lower amount of stellar mass. Note also the strong enhance¬ 
ment of the DM concentration within the central region in 
the SLKA run due to adiabatic contraction (Gnedin et al. 
2004). 

At the near-Eddington growth of the BH in SLDA 
(z = 3.26), the escape velocity profile in the galaxy has 
been boosted by a factor ~ 2, reaching ~ 400 km s -1 at 
r ~ 0.03r v ir (350kms _1 at the galaxy center). These new 
values of the escape velocity in SLDA create a gravita¬ 
tional barrier large enough for the SN winds to be trapped 
within the core of the galaxy, and, therefore, allow mainte¬ 
nance of a steady reservoir of cold gas for the BH to grow 
from. The SLKA escape profile shows a small variation with 
the z = 4.71 profile, it is simply rescaled by a ^ 1 fac¬ 
tor. Although the stellar mass and DM mass content within 
the galaxy bulge are now similar in SLKA and SLDA, the 
amount of gas is still a factor 5 different, and the gas frac¬ 
tion in SLDA is close to 30% while it is below 10% in the 
SLKA. 


3.3 The triggering of nuclear inflows by a major 
merger 

In the simulations with delayed-cooling SN feedback, the BH 
enters its Eddington-limited growth phase once the galaxy 
encounters a significant merger with a stellar mass ratio of 
~1:3. Fig. 8 shows the galaxy effective radius in SHDA as a 
function of its stellar mass for different redshifts. The stel¬ 
lar mass of a galaxy is defined by the total amount of stars 
which belong to the same structure using a galaxy finder al¬ 
gorithm. We employ the AdaptaHOP finder (Aubert et al. 
2004; Tweed et al. 2009) onto the distribution of all stars 
in the simulation. For computing the effective radius, we 
project the stellar mass along the 3 cartesian axes of the 
box, measure the half-mass projected radius of each projec¬ 
tion, and take the geometrical average of the three compo¬ 
nents. In the early stage of the galaxy evolution (z > 4) its 
effective radius stay nearly constant even though the galaxy 
stellar mass increases. Once the 1:3 merger happens from 


redshift 



. 4.71 4.41 4.13 3.88 3.65 3.44 3.26 3.08 2.92 2.77 2.64 2.51 



10 10 10 11 
Mgal (M sun ) 

Figure 8. Evolution with redshift of the effective radius and 
the stellar mass of the central galaxy in SHDA at different red- 
shifts (indicated by the colors). The 1:3 merger is happening from 
z = 3.7 to 3.5 when the effective radius quickly increases and as 
quickly contracts. 



Figure 10. Evolution with redshift of the bar strength of the 
central galaxy in SHDA (black) and SLDA (grey). The bar be¬ 
comes stronger between 3.1 < z < 3.4 after the major merger 
(occurring at 3.4 < z < 3.7) in SHDA, and is the signature of a 
nuclear inflow of gas that triggers the steady Eddington-limited 
growth of the central BH at this redshift. 

z = 3.7toz = 3.5 the galaxy radius quickly increases, but 
then contracts further. 

Fig. 9 illustrates further the connection between the BH 
feeding and the morphological transformation of the galaxy 
in SHDA. This figure shows (at the top) the BH growth rel¬ 
ative to Eddington as a function of redshift with the two 
largest merger mass ratios highlighted in red, and the pan¬ 
els show the stellar emission in rest-frame u, g, i filters at 
various redshifts. We see that above z > 3.6 the galaxy is a 
very disturbed object with no bright galaxy bulge, plenty of 
blue stellar clusters and a disc scale radius almost constant. 
Between 3.21 < z < 3.6, a stellar disc forms with a central 
bar that trigger the inflow of gas in the center of the galaxy 
resulting in a compact bulge later on at z < 3.21. In this 
redshift range, the galaxy settles into a disc-like shape and 
has increased its radius (see also Fig. 8) at the onset of the 
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Figure 9. Top: Eddington ratio of the central BH as a function of redshift for the SHDA. Important mergers for the BH growth are 
indicated as the red horizontal lines with their stellar mass ratio. Small vertical bars pinpoints the redshifts of the stellar images. Bottom: 
Stellar emission of the central galaxy at different redshifts for SHDA as they would be observed face-on through u-g-i filter bands in the 
rest-frame. Extinction by dust is not taken into account. We see three main phases of galaxy evolution: the morphologically disturbed 
phase (z > 3.6), the disc settling phase with a bar (2.77 ^ z ^ 3.6), and the red nugget (z ^ 2.6). 
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merger with its massive companion (that one can see at the 
east/south-east of the panel corresponding to z = 3.6). We 
suggest that the central bar and the resulting nuclear inflow 
is induced by the new gravitational torques exerted by the 
merger, which explain the rise in the BH growth at the near- 
Eddington rate at this redshift. Finally, this merger has al¬ 
lowed sufficient amounts of mass to concentrate around the 
BH through new gas inflows. In return, the gravitational 
barrier (or escape velocity) of the bulge has become large 
enough to prevent the SN winds from depleting the cold 
gas content in the center of the galaxy. The consequence for 
the BH is that there is now plenty of gas residing in the 
bulge, and the BH can, from z — 3.6, grow at the near- 
Eddington rate until it reaches self-regulation z ~ 2.5. At 
the end of the merger (z = 3.4), the galaxy starts to shrink 
and later on (z = 2 . 6 ) form a red nugget with a smooth 
red disc component and a compact stellar bulge. In Fig. B1 
of the Appendix B, we show the morphological evolution 
of the central galaxy in the SLDA run extended to lower 
redshift z — 2. The morphological evolution is remarkably 
comparable to the high-resolution run SHDA, and shows the 
formation of a spiral blue disc galaxy below the redshift of 
z < 2.5. 

In Fig. 10, we show the bar strength of simulations 
SHDA and SLDA. The bar strength is obtained by mea¬ 
suring the inertia tensor of young stars (younger than 10 
Myr) within 1 kpc radius from the center of the galaxy. The 
inertia tensor is then diagonalized to get the three principal 
axis of the ellipsoid e\ < < es- We, then, use those three 

axis to evaluate the strength of the asymmetric distribution 
of young stars (i.e. the so-called “bar”) with the bar strength 
being equal to {e^/e 2 — l)e^/e\. We see that the bar strength 
in SHDA has two main peaks. The first one at z ~ 3.9, which 
correspond to the triggering of accretion bursts onto the cen¬ 
tral BH, though they are discontinuous and not at Edding¬ 
ton (see top panel of Fig. 9 between 3.6 < z < 3.9). The 
second more extended peak at 3.1 < z < 3.4, and that is the 
result of the major merger occurring before at 3.4 < z < 3.7 
and which is responsible for the global destabilization of the 
disc and the central nuclear inflow. This nuclear inflow takes 
the form of a bar that can be seen in the stellar light distri¬ 
bution in Fig. 9. Note that the same bar strength increase is 
observed for the lower resolution run SLDA, except that it 
appears slightly earlier in time (3.1 < z < 3.8) with a lower 
amplitude (that can probably be attributed to resolution ef¬ 
fects). It also explains why the Eddington-limited growth in 
SLDA appears slightly earlier than in SHDA (see top panel 
of Fig. 2). 

Fig. 11 shows the galaxy stellar mass as a function of 
time for the different simulations together with the bulge 
and disc mass content. To compute the bulge and disc mass, 
we build the one-dimensional stellar surface density profile 
extracted from the galaxy finder (in order to remove stellar 
sub-structures using AdaptaHOP), and fit two exponen¬ 
tials £ b ,d = L! 0 ,M exp(-r/rb,d), where E 0 ,M an d n>,d are 
the free parameters of the fit. The bulge Mb and disc Md 
masses are simply Mb,d = Fb,d27rr 2 where the bulge is 
represented by the profile with the smaller rt>. The galaxy 
mass is smaller in all delayed cooling SN feedback runs (SLD, 
SLDA, SHDA) at early times compared to the no feedback 
run (SL), until z ~ 3, but this is essentially due to the sup¬ 
pression of the bulge component. The disc component in SL, 
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Figure 11. Galaxy mass (top panels), bulge mass (middle pan¬ 
els), and disc mass (bottom panels) as a function of time (or 
redshift). In the presence of delayed cooling SN feedback (SLD, 
SLDA and SHDA) the bulge is suppressed at early times (above 
z > 3). The presence of kinetic SN feedback has little impact on 
galaxy mass, both bulge and disc (SLK compared to SL). However 
the presence of the AGN greatly reduces the bulge mass content 
independently both of the SN implementation and resolution. 

SLDA and SHDA is comparable and somewhat larger than 
the disc component of the SL run. 

The kinetic mode of SN feedback, instead, has very lit¬ 
tle impact on the stellar mass content, both for the bulge 
and disc: the mass evolution of the SLK and SL are very 
close. We see that the bulge mass in the SLDA and SHDA 
runs finally settles at z — 3 once the major merger has hap¬ 
pened at z — 3.5. Before the settling of the bulge, its mass 
shows extreme variations with time. This effect is the con¬ 
sequence of the quick removal of gas due to SN explosions, 
creating rapid variations in the bulge gravitational potential 
well and leading to the flattening of the stellar mass distri¬ 
bution (Governato et al. 2010) with the destruction of the 
central stellar clump. Note that the presence of the AGN 
feedback has a non-negligible effect on the mass content of 
the bulge. Once the BH is self-regulating (z ~ 6 for SLKA 
and SHKA, and z ~ 2.5 for SLDA and SHDA), the bulge 
mass growth is regulated by the presence of the AGN feed¬ 
back with more than a factor two difference at z = 2 com¬ 
pared to the simulations without AGN. 

At z = 2, the halo mass is Mh = 10 12 M©, and the 
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time (Gyr) 

Figure 12. BH mass over bulge mass (top panel), BH mass over 
galaxy mass (solid) or disc mass (dotted, bottom panel) as a 
function of time (or redshift) for the SHDA (black), SLDA (grey), 
SHKA (red), and SLKA (orange) simulations. The BH over bulge 
mass ratio is ~ 7 X 10 —4 for all simulations when the BH has 
reached the self-regulation and stays approximately constant. 


expected stellar mass at this redshift obtained from abun¬ 
dance matching is M ga i ~2x 10 10 M© (Moster et al. 2013). 
In all our simulations, the central galaxy ends up being more 
massive than expected from semi-empirical predictions by a 
factor 2 (M ga i ~ 4 x 10 10 M©) for SLKA and by a factor 
3.5 (Mg a i ~ 7 x 10 10 M©) for SLDA. This tension is pos¬ 
sibly due to the missing extra sources of stellar feedback 
including stellar winds and radiation pressure from young 
stars (Hopkins et al. 2011). Including more stellar feedback 
and having a larger star formation efficiency should reduce 
the stellar mass content of the galaxy and delay the forma¬ 
tion of the bulge (Governato et al. 2010; Roskar et al. 2014; 
Agertz & Kravtsov 2015), and, as consequence, further in¬ 
hibits the growth of the central BH. It is worth noting that 
the SN-quenching of BH growth in SLDA occurs until the 
halo reaches a mass of Mh ~ 3 x 10 11 M© at z ~ 3, which 
is smaller than the expected peak of the stellar conversion 
efficiency at z = 3 (Mh, pea k — 3 x 10 12 M©, Moster et al. 
2013). Therefore, if more realistic stellar feedback can pro¬ 
duce lower stellar mass at a corresponding halo mass before 
that peak, we can assume that the halo mass at which this 
SN-quenching of BH growth mechanism becomes ineffective, 
appears for more massive halos, probably around the peak 
of stellar conversion efficiency. 

In the case with delayed cooling SN feedback, the BH 
starts to grow (see Fig. 2) when the bulge definitively set¬ 


tles at z ~ 3.5, while in the kinetic feedback case bulge 
and BH start growing earlier on. In Fig. 12 we compare, in 
the top panel, the ratio between BH and bulge mass as a 
function of time. In SLKA and SHKA a roughly constant 
ratio, ~ 7 x 10 - , is maintained from early times, z ~ 6 
all the way to the end of the simulations. This value is only 
slightly lower than the ratio found at z = 0, 1 — 3 x 10 -3 . 
SHDA and SLDA present instead a very different evolu¬ 
tion: the BH’s growth is very slow until z ~ 3.5, and the 
bulge has not settled yet. Therefore, before that time, the 
ratio between BH and bulge mass varies rapidly, despite the 
BH mass remaining roughly constant. In this case a mea¬ 
sure of the BH-to-bulge ratio would indicate a value typi¬ 
cally higher than ~ 10 -3 . The bulge is able to grow, and to 
let the BH grow, only after the potential of the bulge has 
reached a large enough value (cf. section 3.2). We have also 
evaluated the ratio of the BH to total galaxy mass (bot¬ 
tom panel of Fig. 2). In fact, except for a few exceptions 
(Peng et al. 2006; Schramm & Silverman 2013) high-z mea¬ 
surements of the masses of quasar and AGN hosts cannot 
perform a robust bulge-to-disc decomposition, and often the 
total stellar mass is measured (but compared to the BH- 
bulge relationship at z = 0). In SLKA and SHKA the ratio 
gently decreases with cosmic time, highlighting that the to¬ 
tal stellar mass grows faster than the BH. Therefore, in the 
SLKA and SHKA cases, the BH bulge-mass relation would 
be roughly independent of redshift below z < 6 while the 
BH-galaxy mass relation would decrease smoothly with red¬ 
shift. In SHDA and SLDA, the situation is more complex: 
initially the galaxy grows faster than the BH (z > 3.5), then 
the BH grows faster than the galaxy (3.5 < z < 2.5), and 
eventually they start tracking each other. In this case, low 
mass BHs at high-redshift would be considered overmassive 
with respect to their bulge, but undermassive if the mea¬ 
sured quantity is the total galaxy mass, and the comparison 
is made to today’s BH-bulge relationship (i.e., an expecta¬ 
tion value of ~ 10 -3 ). 

3.4 Extension to the surrounding galaxies 

To check whether the impact of efficient SN explosions ap¬ 
plies to low-mass galaxies in a general fashion, we measure 
the BH mass and galaxy mass relation for the other galaxies 
within the zoom region. We only consider galaxies that are 
at a distance smaller than 5r v ir of the main galaxy. Fig. 13 
shows the values of MBH-M ga i for these galaxies at various 
redshifts in the SLDA and SLKA runs, with the observed 
MBH-M ga i relation at z — 0 from Haring & Rix (2004). As 
already highlighted in the previous sections for the main 
galaxy, we recognize the slowing down of the growth of the 
most massive BH (the BH of the main galaxy) at high red¬ 
shift in the run with delayed cooling SN feedback (SLDA), 
while the most massive BH grows quickly together with 
its host galaxy mass in the kinetic SN run (SLKA). Other 
smaller mass galaxies behave in a similar fashion: in SLKA 
the BH mass follows the MBH-M ga i relation (though with 
a larger normalization factor) and BHs starts growing in 
galaxies of M ga i ~ 10 8 M©, in SLDA the BH growth is de¬ 
layed in low-mass galaxies, i.e. BHs start growing in more 
massive galaxies of M ga i ~ 2 x 10 9 M©, the normalization 
of MBH-M ga i relation is further increased. Note that there 
is a plateau in our MBH-M ga i relation at low galaxy masses 
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Figure 13. BH mass versus stellar mass at different redshifts as indicated in the different panels for SLDA (grey circles) and SLKA (orange 
stars). The main central galaxy is indicated with the filled symbol. The solid line is for the z = 0 observed relation from Kormendy & Ho 
(2013). In the SLKA simulation, BHs reach the self-regulated state very early-on, while the presence of delayed-cooling SN feedback in 
SLDA prevents the BHs in the less gravitationally bound bulges from growing rapidly. 


which is due to our choice of initial seed BH mass 10 5 Mq 
( there are BHs with smaller values because we cannot take 
more mass within a gas cell that is actually available), there¬ 
fore, those BHs are already overly massive for their galaxy 
mass. 

In Fig. 14, we compare the projected gas densities of 
a few galaxies in SLDA and SLKA at 2 = 2. Galaxies in 
SLDA exhibit a porous and turbulent ISM that is driven by 
SN explosions. In comparison, galaxies in SLKA systemati¬ 
cally show a more uniform gas distribution with spiral arms 
encompassing star-forming gas clumps but the ISM is not 
as porous as in the SLDA run. As discussed earlier for the 
main galaxy, this extra turbulence (or porosity) driven by ef¬ 
ficient SN explosions, is responsible for the inhibition of BH 
growth in the low-mass galaxies. Therefore, it seems plausi¬ 
ble that this inhibition mechanism should hold on average 
for all galaxies with galaxy masses below a few 10 9 M©. 


4 CONCLUSIONS 

Our set of cosmological hydrodynamical simulations (the 
Seth suit of simulations) including either SN feedback with 
a kinetic model from Dubois & Teyssier (2008) or a delayed 
cooling prescription from Teyssier et al. (2013) shows that 
the strength of SN feedback can dramatically alter the BH 
evolution at early epochs of galaxy formation, or conversely 
in low-mass galaxies. In the former (weak SN feedback) case, 
a cold gas reservoir is maintained in the center of the galaxy 
sustained by inflows of dense star-forming clumps of gas, 
which provides an efficient gas fuelling of the central BH, 
driving an Eddington-limited growth early on. In the later 
(strong SN feedback) case, the energy release from SNe is 


maximised, which allows for the destruction of the dense 
gas clumps within the core of the galaxy. In return, the 
early growth of the BH is prevented due to the lack of a 
sustained cold gas reservoir, and this effect is systematically 
observed for all galaxies within the simulated region. We 
have shown that the characteristic SN wind velocity is of 
the order of 270kms _1 , a value larger than the escape ve¬ 
locity created by the gravitational potential of a galactic 
bulge of Mb < 10 9 M© (with a fixed radius of 100 pc). This 
efficient mode of SN feedback prevents the rapid growth of 
the bulge although the disc component continues growing al¬ 
most unimpeded. Once enough mass has accumulated in the 
center of the galaxy, in our case driven by a major merger, 
the gravitational barrier prevents the removal of the cold 
gas and the BH can grow at a near-Eddington rate until it 
reaches self-regulation. 

For the efficient SN feedback scenario, the bulge mass 
~ a few 10 9 M© above which the accretion at the Edding¬ 
ton rate turns on is a factor ~ 10 smaller than the bulge 
mass ~ 10 10 M© below which abnormally light BHs are ob¬ 
served and exhibit a larger scatter in the BH-bulge correla¬ 
tion (Greene et al. 2010; Jiang et al. 2011; Graham & Scott 
2015). In our simulations, only feedback from type II SNe 
has been included. An extra source of energy or momentum 
input such has radiation pressure from young stars could en¬ 
hance the wind velocities, and, thus, increase the bulge mass 
above which the gas ejected by stellar feedback is trapped in 
the bulge. Therefore, the fact that BHs undergo less growth 
in low-mass galaxies could be the signature of a strong stel¬ 
lar feedback activity, which is able to efficiently destroy the 
cold clumps of star-forming gas, establishing an indirect ev¬ 
idence of short-lived clusters of gas below a given mass. 

We have tested those two scenarios with respect to spa- 
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Figure 14. Projected gas densities for three other galaxies than 
the central one at z = 2 for SLDA (left column) and SLKA (right 
row). The images are 10 kpc across. These galaxies are not sub¬ 
structures of the main halo, and each of this image in SLDA cor¬ 
respond to its equivalent in SLKA (i.e. same halo host). Galaxies 
simulated with delayed cooling SN feedback appear more dis¬ 
turbed and turbulent than galaxies simulated with kinetic SN 
feedback, explaining the lack of early BH growth that is a generic 
effect in low-mass galaxies. 


tial resolution effects at Ax = 34.8 pc and Ax — 8.7, and 
our results are qualitatively independent on the adopted res¬ 
olution. The Eddington growth phases appear at the same 
epochs: as soon as the galaxy forms for the case of kinetic SN 
feedback runs, and once the major merger happens in the 
delayed-cooling SN feedback model. Even the quantitative 
aspects show relatively good agreement: the bulge and disc 
mass are in remarkable good agreement, and the BH mass 
at self-regulation end up being the same within a factor 2. 

In the Seth simulation, the 1:3 major gas-rich merger 
happening at z — 3.5 is the culprit for the triggering of the 
BH growth in the runs including the delayed cooling pre¬ 
scription for SN feedback. It does not imply that all galaxies 
would necessarily start growing their BH every time they en¬ 
dure a major merger, because they probably need to reach 
a mass sufficiently large enough (Mb > 10 9 M©, therefore 


Mg a i > 10 9 M©). Once the galaxy is massive enough, a per¬ 
sistent bulge needs to be built, and the major merger is 
the leading mechanism for the Seth simulation because of 
the induced torques and consequent nuclear inflows. Though 
major mergers are the favourite candidates for creating the 
strongest torques, one can imagine that a massive flow of 
cold streams misaligned with the galactic spin would have 
a similar impact on the galaxy. Only extra simulations with 
different accretion histories, and in particular more quies¬ 
cent objects, could help us to fully address that question and 
verify whether major mergers are systematically responsible 
for the ignition of BH growth in the case of efficient SN 
feedback. 

The detailed interaction between star formation, SN 
feedback, BH growth and AGN feedback is probably the 
key to understanding the establishment of the scaling re¬ 
lation between BHs and their host galaxies, their redshift 
evolution, the overall co-evolution of BHs and galaxies, and 
as a consequence the properties of the galaxy population. 
At high redshift, estimates of BH and galaxy stellar masses 
of AGN samples typically find that the BH mass is larger 
than expected from this correlation at a given stellar mass 4 , 
except for submillimeter-selected galaxy (Alexander et al. 
2008). Willott et al. (2010), however, based on the compar¬ 
ison between the BH and galaxy mass functions at z ~ 6, 
suggest that most high-z galaxies host BHs that are smaller 
than expected from this correlation at a given stellar mass 
(or host no BH at all). Similar conclusions are reached with 
X-ray searches of z > 6 galaxies: current data are sensitive 
to 10 5 —10 6 M© BHs accreting at 1 — 10% Eddington (Willott 
2011; Fiore et al. 2012; Cowie et al. 2012; Treister et al. 
2013; Giallongo et al. 2015; Weigel et al. 2015), but very 
few, if any are found, implying strong limits on the early 
growth of BHs in relatively small galaxies. The observational 
evidence, therefore, seems to suggest that high-redshift BHs 
in the most massive galaxies are actively accreting and are 
on or above the correlation with their hosts, while at lower 
galaxy mass the behaviour is opposite (Volonteri & Stark 
2011). This new SN-regulated BH growth mechanism is a 
promising tool to explain why BHs in relatively low-mass 
galaxies at high-redshift are unable to grow. 

Though our method for modelling SN feedback is still 
very crude, it highlights two extreme regimes. Future im¬ 
provements should model SNe with more physically moti¬ 
vated assumptions using the correct amount of momentum 
imparted by the explosions (Agertz et al. 2013; Kimm et al. 
2015) and study the interplay between circumnuclear star 
formation and black hole growth. We should also consider 
the steady input of energy from stellar winds and radia¬ 
tion from bright OB associations (e.g. Hopkins et al. 2011) 
with dedicated simulations to probe their effect on the BH 
growth. Note that, even though there will be extra energy 
delivered by these new stellar feedback processes, it is un¬ 
clear how they will interact with the gas content and star 
formation. As a matter of fact, stellar winds and UV radia¬ 
tion from young stars initiate earlier than type II SNe, and 
could result in an earlier impact on the gas content of the 


4 High-redshift stellar masses are not equivalent to bulge masses, 
nor are measured in a self-consistent way throughout different 
samples and at different redshift. 
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star-forming cloud but could also limit the total amount of 
stars formed in one star-formation cycle, thus, limiting the 
impact of feedback (Rosdahl et al. 2015). However, a sce¬ 
nario with efficient SN (or stellar) feedback mechanism is 
probably favoured over a less extreme mode in order to ex¬ 
plain the low stellar mass content (e.g. Hopkins et al. 2014) 
and the absence of bulge in dwarf galaxies (Governato et al. 
2010). We dedicate these investigations to future work. 

We have limited our simulations to the high-redshift 
(z > 2 ) environment of the progenitor of a group of galaxies 
and have shown that all low-mass galaxies (M ga i < 10 9 M©) 
within this region exhibit un-grown BHs. Though, we have 
not explicitly demonstrated that the same mechanism is at 
play in the low-redshift dwarf galaxies, we argue that they 
should behave similarly, particularly given that galaxies at 
low redshift are less compact than their high-redshift equiva¬ 
lent, and, thus, it is easier for SN explosions to drive outflows 
from them. This epoch of galaxy evolution remains to be in¬ 
vestigated and our proposed mechanism to be validated by 
new numerical experiments of low-redshift dwarf galaxies. 
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APPENDIX A: SETTING UP PARAMETERS 
FOR THE DELAYED COOLING MODE OF SN 
FEEDBACK 


The time variation of the non-thermal energy component in 
star-forming regions follows 


DeNT . 10 51 erg cnt 

Di = VsNPs 10M© " 


(Al) 


with the first term of the right-hand side being the injec¬ 
tion of energy by SN explosions, and the second term is the 
decay rate due to the dissipation of the energy by sub-grid 


processes (decay of turbulence, cooling of cosmic rays, etc.). 
At equilibrium the two terms exactly compensate and lead 
to 


eNT = rjSNps 


10 51 erg 
10 M© 


tdh 


(A2) 


ONT 

p 


■ ?7SN- 


10 51 erg Idiss 

10M 0 tB(p) f 


(A3) 


for a Schmidt law. It comes that the corresponding non- 
thermal velocity dispersion is (eNT = 0.5pcrN T ): 


O' NT — 2?7SN 


10 51 erg l diss 

10 M 0 


(A4) 


The injected energy from SNe is transferred from unresolved 
small scales, with high gas density where it is deposited, up 
to resolved larger scales, with low gas density. Explicit re¬ 
solved turbulence will dissipate over its own time-scale which 
is provided by the dynamics of shocks in the ISM. Note that, 
here, we derive a model for sub-grid unresolved turbulence, 
and the time-scale to dissipate this energy component (non- 
thermal energy) is turbulent crossing-time over a few res¬ 
olution elements (typically one Jeans length) lj — N ce \\Ax 
(with iV C eii ^ 4), so that tdiss — ^j/ctnt- Therefore, eq. (A4) 
is replaced by 


ctnt 



10 51 erg lj \ 
10 M 0 ts{p) e *) 


1/3 


(AS) 


We can simplify the expression further by imposing that the 
gas free-fall time is that at the threshold of star formation. 
Thus 


t*(p) = 


ts(p) 

e* 


1 / 3tt 

Z V 32Gpo ' 


(A6) 


Finally, we can write ctnt as a simple function of all known 
parameters for star formation and SN feedback efficiency, 
and resolution: 



where no is the number density threshold for star formation 
(no = nu,o/Xu, Xu = 0.76). The dissipative time-scale is 
also fully determined by the choice of p sn, e*, Ax and no: 



For this set of simulations we have taken V ce ii = 4 corre¬ 
sponding to the typical Jeans length at high gas densities 
(pressure-supported by the polytropic equation of state, see 
section 2.2). Note that Xnt becomes resolution independent 
for our given choice of polytropic equation of state. Since 
4 x Ax — lj oc y/To /no, the density dependancy in equa¬ 
tion (A7) vanishes. This simple model ensures that SN en¬ 
ergy deposit can support the propagation of the blast wave 
up to iV C eiiAx without any radiative loss since cooling is 
blocked up to that distance. Beyond that distance, we con¬ 
sider that the blast wave is resolved and radiative losses can 
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proceed, which is ensured by the calibration of the dissipa¬ 
tion time-scale tdiss to that of the sound-crossing time up to 
N ce \\Ax for an energy deposit of 0.5/9 <Tn T . 


APPENDIX B: GALAXY MORPHOLOGY AT 
LOWER RESOLUTION 

In Fig. Bl, we show the stellar emission in the SLDA run 
with lower resolution (Ax = 34.8 pc) than the SHDA run 
(Ax = 8.7 pc) presented in Fig. 9. The morphology of the 
galaxy follows qualitatively the same behavior than the high- 
resolution run with the morphologically disturbed shape at 
z > 3.6, the disc-like shape with a bar component (see e.g. 
panels z = 3.44 and 3.26 in Fig. Bl) during 3.1 < z < 3.8 
(see Fig. 10), the red nugget phase at 2.53 ^ z ^ 2.77, and 
the formation of a spiral disc galaxy at z ^ 2.45 with a 
red bulge. As for the SHDA run, the near-Eddington BH 
growth starts around z = 3.5, once the galaxy quits the 
morphological disturbed phase and that the massive merger 
triggers the formation of a bar and nuclear inflows. 
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Figure Bl. Same as Fig. 9 for the low-resolution run (SLDA). We have extended the stellar images to lower redshifts that we have been 
able to reach due to the lower resolution of the simulation. Note that though the resolution is now Ax = 34.8 pc, the galaxy looks alike 
the SHDA run (Fig. 9). We also see the formation of a disc spiral galaxy in the redshift range not probed by SHDA (z ^ 2.31). 
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